///////////////////////////////////////////////////////////////////////////////////////////////////
/// \file vmc_walker.h
///
/// \brief   header file for the variational monte carlo walker class
///
///////////////////////////////////////////////////////////////////////////////////////////////////

#ifndef PCPS_VMC_WALKER_HEADER
#define PCPS_VMC_WALKER_HEADER

#include<src/pcps.h>
#include<src/input.h>
#include<src/walker.h>
#include<src/configuration.h>
#include<src/accumulator.h>
#include<src/coeff_obj.h>
#include<src/restricted_pairing_function.h>

namespace pcps {

  ///////////////////////////////////////////////////////////////////////////////////////////////////
  /// \brief  A class for a variational monte carlo walker
  ///
  ///////////////////////////////////////////////////////////////////////////////////////////////////
  template<class S> class VMC_Walker : public pcps::Walker {

    private:

      /// \brief  whether the most recently proposed configuration was accepted
      bool _accepted;

      /// \brief  a counter to track how long the walker has been on the current configuration
      double _count;

      /// \brief  the number of steps taken since the last accumulation
      int _steps;

      /// \brief  the configuration the walker is on
      pcps::Configuration<S> _config;

      /// \brief  the proposed configuration to move to next
      pcps::Configuration<S> _next_config;

      /// \brief  an object specifying the difference between the current and proposed configurations
      pcps::ConfigDiff<S> _config_diff;

      /// \brief  a random number generator
      pcps::LaggedFibonacci<> _lfg;

      /// \brief  pointer to the coefficient object used as the guiding function
      boost::shared_ptr<pcps::CoeffObj<S> > _gf;

      /// \brief  pointers to coefficient objects that should be updated as we move around
      std::vector<boost::shared_ptr<pcps::CoeffObj<S> > > _co_vec;

      /// \brief  pointers to accumulators
      std::vector<boost::shared_ptr<pcps::Accumulator<S> > > _accumulators;

    public:

      ///////////////////////////////////////////////////////////////////////////////////////////////////
      /// \brief Constructs the walker.
      ///
      /// \param[in]      userinp  the user's input options
      ///
      ///////////////////////////////////////////////////////////////////////////////////////////////////
      VMC_Walker(const pcps::Input & userinp)
        : _accepted(false),
          _count(1.0),
          _steps(0),
          _config(userinp),
          _next_config(userinp),
          _config_diff(userinp),
          _lfg(),
          _gf(),
          _co_vec(),
          _accumulators()
      {}

      /// \brief  returns a reference to the walker's configuration
      pcps::Configuration<S> & config() { return this->_config; }

      /// \brief  returns a constant reference to the walker's configuration
      const pcps::Configuration<S> & config() const { return this->_config; }

      /// \brief  returns a reference to the walker's random number generator
      pcps::LaggedFibonacci<> & lfg() { return this->_lfg; }

      /// \brief  returns a constant reference to the walker's random number generator
      const pcps::LaggedFibonacci<> & lfg() const { return this->_lfg; }

      /// \brief  returns a reference to the walker's guiding function pointer
      boost::shared_ptr<pcps::CoeffObj<S> > & gf() { return this->_gf; }

      /// \brief  returns a constant reference to the walker's guiding function pointer
      const boost::shared_ptr<pcps::CoeffObj<S> > & gf() const { return this->_gf; }

      /// \brief  returns a reference to the walker's vector of other (non-guiding) coefficient objects
      std::vector<boost::shared_ptr<pcps::CoeffObj<S> > > & co_vec() { return this->_co_vec; }

      /// \brief  returns a constant reference to the walker's vector of other (non-guiding) coefficient objects
      const std::vector<boost::shared_ptr<pcps::CoeffObj<S> > > & co_vec() const { return this->_co_vec; }

      /// \brief  returns a reference to the walker's vector of accumulators
      std::vector<boost::shared_ptr<pcps::Accumulator<S> > > & accumulators() { return this->_accumulators; }

      /// \brief  returns a constant reference to the walker's vector of accumulators
      const std::vector<boost::shared_ptr<pcps::Accumulator<S> > > & accumulators() const { return this->_accumulators; }

      ///////////////////////////////////////////////////////////////////////////////////////////////////
      /// \brief Proposes a move to a new configuration, which is assumed to be stored internally.
      ///
      /// \param[in]      userinp  the user's input options
      ///
      ///////////////////////////////////////////////////////////////////////////////////////////////////
      void propose_move(const pcps::Input & userinp) {

        // propose where to move to and get the proposal density ratio
        const double pdr = pcps::Configuration<S>::propose_move(userinp, _lfg, _config, _next_config, _config_diff);

        // calculate the ratio of the guiding function's weight for the new and old configurations
        S gf_ratio = ( _gf->get_value(_next_config, _config, _config_diff) / _gf->get_value() ).get();

        //const S gf_new_val = gf_ratio * _gf->get_value().get();
        //const S gf_old_val = _gf->get_value().get();

        const pcps::RestPairingFunc<S> * const rpf = dynamic_cast<pcps::RestPairingFunc<S> *>( _gf.get() );
        if ( !rpf->wfn_sampling() )
          gf_ratio = pcps::unity<S>() * std::sqrt( std::sqrt( pcps::square_norm(gf_ratio) ) );

        // choose a random number and decide whether to accept the proposed move according to the Metropolis algorithm
        _accepted = ( pcps::uni_01(_lfg) < pdr * pcps::square_norm(gf_ratio) );

//        // check the value of the old configuration
//        if ( false && !rpf->wfn_sampling() ) {
//          const S old_val = _gf->get_value().get();
//          const S cor_val = rpf->correct_value(_config, 'v');
//          //std::cout << boost::format("    proposed  %20.12f          correct  %20.12f           %s")
//          //             % old_val
//          //             % cor_val
//          //             % ( _accepted ? "accepted" : "" )
//          //          << std::endl;
//          assert( std::abs( old_val - cor_val ) < 1.0e-9 );
//        }
//
//        // check the value of the new configuration
//        if ( false && !rpf->wfn_sampling() && _accepted ) {
//          const S cor_val = rpf->correct_value(_next_config, 'v');
//          std::cout << boost::format("  proposed  %20.12f     correct  %20.12f           %s")
//                       % gf_new_val
//                       % cor_val
//                       % ( _accepted ? "accepted" : "" )
//                    ;
//                    //<< std::endl;
//          assert( std::abs( gf_new_val - cor_val ) < 1.0e-9 );
//        }
//
//        if (true) {
//          for (int i = 0; i < _config.fock_config().size(); i+=2)
//            std::cout << _config.fock_config().at(i);
//          std::cout << " ";
//          for (int i = 1; i < _config.fock_config().size(); i+=2)
//            std::cout << _config.fock_config().at(i);
//          std::cout << boost::format("    proposed  %20.12f          current  %20.12f           %s")
//                       //% (dynamic_cast<RestPairingFunc<S> *>(_gf.get()))->prop_ratio_sum()
//                       % gf_new_val
//                       % gf_old_val
//                       % ( _accepted ? "accepted" : "" )
//                    << std::endl;
//        }

      }

      ///////////////////////////////////////////////////////////////////////////////////////////////////
      /// \brief Accumulate data at the current configuration.
      ///
      /// \param[in]      userinp  the user's input options
      ///
      ///////////////////////////////////////////////////////////////////////////////////////////////////
      void accumulate(const pcps::Input & userinp) {

        //if ( _steps != 0 )
        //  return;

        // if the move was accepted, accumulate data for the configuration the walker is about to leave
        if (_accepted) {
          typename std::vector<boost::shared_ptr<pcps::Accumulator<S> > >::iterator it = _accumulators.begin();
          for ( ; it != _accumulators.end(); it++)
            (**it).accumulate(userinp, _count, _config);
        }

      }

      ///////////////////////////////////////////////////////////////////////////////////////////////////
      /// \brief Move the walker according to the proposed move and whether it was accepted.
      ///
      /// \param[in]      userinp  the user's input options
      ///
      ///////////////////////////////////////////////////////////////////////////////////////////////////
      void move(const pcps::Input & userinp) {

        // if the move was accepted...
        if (_accepted) {

          // update the guiding function for the new configuration
          _gf->update(_next_config, _config, _config_diff);

          //std::cout << boost::format("move accepted,  new coeff = %20.12f        noiseless coeff = %20.12f")
          //             % _gf->get_value().get()
          //             % _gf->noiseless_value(_next_config).get()
          //          << std::endl
          //          << std::endl;

          // update the non-guiding coefficient value objects for the new configuration
          typename std::vector<boost::shared_ptr<pcps::CoeffObj<S> > >::iterator it = _co_vec.begin();
          for ( ; it != _co_vec.end(); it++)
            (**it).update(_next_config, _config, _config_diff);

          // update the current configuration to reflect the move
          _config = _next_config;

          // reset the counter for how long we have been on the current configuration
          _count = 1.0;

        // if the move was rejected...
        } else {

          //std::cout << "move rejected" << std::endl << std::endl;

          // update the counter for how long we have been on the current configuration
          _count += 1.0;

        }

        //_steps = ( _steps + 1 ) % userinp.nparticles();
        //_count = 1.0;

      }

      ///////////////////////////////////////////////////////////////////////////////////////////////////
      /// \brief Perform any remaining accumulations and communicate the accumulators' results.
      ///
      /// \param[in]      userinp  the user's input options
      ///
      ///////////////////////////////////////////////////////////////////////////////////////////////////
      void final_accumulate(const pcps::Input & userinp) {

        // accumulate data for the final configuration
        for (typename std::vector<boost::shared_ptr<pcps::Accumulator<S> > >::iterator it = _accumulators.begin();
             it != _accumulators.end();
             it++)
          (**it).accumulate(userinp, _count - 1.0, _config);
          //(**it).accumulate(userinp, _count, _config);

        // finilize the accumulators
        for (typename std::vector<boost::shared_ptr<pcps::Accumulator<S> > >::iterator it = _accumulators.begin();
             it != _accumulators.end();
             it++)
          (**it).finalize(userinp);

        //std::cout << std::endl;

      }

      ///////////////////////////////////////////////////////////////////////////////////////////////////
      /// \brief reset the walker in preparation for sampling
      ///
      /// \param[in]      userinp  the user's input options
      ///
      ///////////////////////////////////////////////////////////////////////////////////////////////////
      void reset(const pcps::Input & userinp, const bool reset_co = true) {

        if (reset_co) {

          // initialize the guiding function
          _gf->initialize(_config);

          // initialize the non-guiding coefficient value objects
          for (typename std::vector<boost::shared_ptr<pcps::CoeffObj<S> > >::iterator it = _co_vec.begin();
               it != _co_vec.end();
               it++)
            (**it).initialize(_config);

        }

        // reset the accumulators
        for (typename std::vector<boost::shared_ptr<pcps::Accumulator<S> > >::iterator it = _accumulators.begin();
             it != _accumulators.end();
             it++)
          (**it).reset(_config);

        // initialize the count of how long we have been at the current configuration
        _count = 1.0;

        // reset the number of steps taken since the last accumulation
        _steps = 0;

      }

  };

}

#endif
